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The interaction and exchange-correlation contributions to the ground-state energy of an arbitrary 
many-electron system can be obtained from a spherical average of the wavevector-dependent diagonal 
structure factor (SF). We model the continuous-k spherically averaged SF using quantum Monte 
Carlo calculations in finite simulation cells. We thus derive a method that allows to substantially 
reduce the troublesome Coulomb finite-size errors that are usually present in ground-state energy 
calculations. To demonstrate this, we perform variational Monte Carlo calculations of the interaction 
energy of the homogeneous electron gas. The method is, however, equally applicable to arbitrary 
inhomogeneous systems. 

PACS numbers: 71.15.-m,71.10.Ca,71.45.Gm 



X3 

O 

o 



> 
o 

o 
o 



o 



X 
S3 



Together with density- functional theory (DFT) P, 
Quantum Monte Carlo (QMC) calculations belong to the 
bedrock of computational solid state physics One 
major problem, the focus of this paper, encountered in 
QMC solid-state applications are the Coulomb finite-size 
effects. These originate in the periodic Ewald interaction 
that is typically used to model the electron-electron (e-e) 
Coulomb interaction in a periodic geometry When 
using the Ewald interaction, the long-range nature of the 
Coulomb interaction yields spurious contributions to the 
interaction energy caused by the interaction of an elec- 
tron with the periodically repeated copies of its exchange- 
correlation (xc) hole. Such finite-size effects are usually 
dealt with by increasing the system size and monitor- 
ing the convergence of the relevant data. However, the 
Coulomb finite-size errors of the interaction energy are 
known to scale as N being the number of electrons 
in the supercell, and convergence is therefore slow 
An alternative involves replacing the periodic Ewald in- 
teraction by a "model periodic Coulomb" (MPC) inter- 
action that converges faster 0, Q ■ 
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In this Letter, we present an approach that reduces 
finite-size errors by keeping the true Coulomb interaction 
and going to the core of the issue by using QMC to model 
the spherical average of the diagonal structure factor (SF 
of extended systems. This is in contrast to a paper 
we recently learnt of that uses the (non-spherically aver- 
aged) SF to correct QMC data at long wavelengths. Our 
starting point is the link between the interaction energy 
and the spherically averaged wavevector-dependent diag- 
onal SF Sk of an arbitrary many-electron system. Since 
large k wavevectors are related to small electron separa- 
tions, the large- /c behavior of QMC simulations ought to 
be correct. This is not the case for the small-fc behavior 
(long wavelengths), which in the case of a finite simula- 
tion cell yields spurious Ewald interactions. However, the 
long-range (small-fc) behaviour can be obtained differ- 
ently, either from the known constraint Sk~^o = or 



from alternative calculations based on, e.g. the random- 
phase approximation (RPA) which is known to describe 
collective excitations correctly and becomes exact at long 
wavelengths ■ We aim to have the best of both worlds! 

The method we propose bears some resemblance to the 
MPC interaction. However, while the MPC modifies the 
interaction, here we keep the standard Coulomb interac- 
tion and model the QMC correlations in a fc dependent 
way. The advantage is that our method ought to yield 
improvements where others fail e.g. the exchange 
hole is known to be long-ranged, decaying as In 
a finite simulation cell this results in a finite size error 
of the exchange energy cx A^~^/^. The error at large r, 
however, corresponds to Sk at small fc where it can be 
replaced easily by the correct asymptotic value. 

The interaction energy of an arbitrary many-electron 
system is usually expressed as the sum of the Hartree 
energy (which in the case of an infinite Jellium model 
is exactly cancelled by the Coulomb energy due to the 
presence of the positive background) and the so-called xc 
interaction energy Uxc- Uxc corresponds to the attractive 
interaction between each electron and its own xc hole. 
Starting from the spherical average n2,c(r,w) of the xc 
hole density na,c(r, r') at r' around an electron at r, one 
finds El 



N 



dk{Sk-l), 



(1) 



where Sk is the spherical average of the diagonal struc- 
ture factor 5*1; k: 
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n^c(r,it), (2) 



n(r) being the electron density at r. 

Implicitly we are considering an infinite system: fc is 
a continuous variable. However, Sk for a QMC system 
contains irregularly spaced delta peaks that on integra- 
tion give the QMC Uxc- Our model (see Fig. ^) assumes 
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FIG. 1: An electron at « = inside a box of size L is sur- 
rounded by a periodic xc-hole p^c ''^ ■ Uxc is the Coulomb en- 
ergy between the electron and its xc-hole pxc- Our method ef- 
fectively replaces p—p^^'^ (EWALD in the figure) by the den- 
sity p (SF) beyond a given cutoff uq. Clearly, the interaction 
between an electron and its xc hole at u = L, 2L, ... is more 
accurate (closer to REAL) and converges much faster than in 
the standard QMC case (EWALD). Close to the boundary of 
the box pxc and p^c '''^ diverge as both have to integrate to 
— 1. Since p^c in the SF case jumps to at uq this sum-rule 
is violated, but the resulting error can be corrected for easily 
(see text). 



that correlations are non-periodic and beyond a cutoff 
radius uq are due only to variations in the density, so 
that beyond Uq there is no contribution to Uxc- For uq, 
we choose the Wigner-Seitz radius uws of the simulation 
cell, and the structure factor Sk can then be sampled di- 
rectly during the QMC run. We find: 



Sk — sl + sl 
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using the dimensionless quantities fiq — Vuq, k — kuQ, 
and q = \q\uo. f = 47ruo/(3T^) is the volume fraction of 
the super cell that contributes to Sl, V is its total vol- 
ume, and Uq denotes the Fourier transform of the electron 
density n(r). Apart from n(r) the sampling of Sk only 
needs \rj — rj| which is readily available in most QMC 
codes making the SF easy to implement. Sl^ , which is 
due solely to variations in the density, cancels the Hartree 
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FIG. 2: This figure shows the raw QMC and the adjusted 
Hartree-Fock SF for the smallest and the largest system (2*27 
and 2*307 electrons respectively) as well as the exact result. 



contribution to Eq. iQ. Note that the k is continu- 
ous even in the case of a periodic system for which the q 
vectors are discrete. This is deliberate, as the sampling 
of a periodic QMC system models an extended (non pe- 
riodic, continuous k) system. 

Equations Q and ^ do not include the entire xc hole 
in the QMC sampling, as / < 1. As a result, our raw 
SF differs from zero at A: = 0. This corresponds to the 
amount of the xc hole that is missed and which is lo- 
cated in the corners of the simulation cell beyond the 
cutoff radius J5]. However, due to the periodic bound- 
ary conditions the QMC description in these corners is 
unlikely to be accurate, so not much information is lost. 
Below we show that some of the residual error can be 
corrected for easily and efficiently. 

We performed Hartree-Fock (HF) and Variational 
Monte Carlo (VMC) calculations for the homogeneous 
electron gas, using the CASINO package [3. The calcu- 
lations employ plane-wave Slater determinants with and 
without a Jastrow factor. The latter corresponds to HF 
calculations where the exact result is known The 
systems we studied are non-polarised in a face-centered 
cubic simulation cell with the number of electrons rang- 
ing from 2 * 27 to 2 * 307 at Ts = 1. This corresponds to a 
Wigner-Seitz radius uws ranging from 3.420 to 7.689 in 
atomic units, which we use throughout. The interaction 
energy is evaluated using either our SF-based approach 
or the standard Ewald interaction. The SF is sampled at 
1000 equally spaced points ranging from to 10. In the 
case of VMC calculations, the Jastrow factor was con- 
verged using several iterations of variance minimisation. 

Let us first look at a pure Slater determinant of plane 
waves. Fig. |21 plots our HF calculation of the SF next 
to the exact Hartree-Fock SF. Fig. El shows the conver- 
gence of Uxc to the known HF (exchange) energy. The 
Ewald data shows the familiar finite-size errors, while the 
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FIG. 3: Convergence of the potential energy with system 
size (2*27 to 2*307 electrons) using a Slater determinant only. 
The solid line is the HF result for an infinite system. 



FIG. 5: Convergence of the potential energy with system 
size (2*27 to 2*307 electrons) using a Slater-Jastrow wave- 
function. The standard Ewald result is shown next to the 
raw and corrected SF result. 
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FIG. 4: The figure shows the raw QMC and the adjusted 
Slater-Jastrow SF for the smallest and the largest system 
(2*27 and 2*307 electrons respectively). The RPA SF is also 
shown. It seems incorrect for k > 0.5. Hence, only for systems 
large enough such that 0.5 ^ kc would interpolation with the 
RPA SF improve U^c- 

SF yields an interaction (exchange) energy that is essen- 
tially flat, consistent with the elimination of the Coulomb 
finite-size error. 

A similar analysis can be performed with the correlated 
Slater-Jastrow many-electron wavefunction (Figs. ^ and 
O . As in the case of the HF calculation, the Ewald data 
exhibits finite-size errors, which should scale as 1/iV. In 
contrast, our SF-based calculations exhibit a systematic 
error with the opposite sign. The Coulomb finite-size 
error has been eliminated, but we are missing a bit of 
the xc hole (located at the corners of the simulation cell), 
which yields a structure factor Sk that differs from zero 



at fc = resulting in an erroneous xc interaction energy. 

We now analyse the behavior of the structure factor 
Sk at small k. As the system size increases, one would 
expect the SF to improve at small k. Fig. |21 shows our 
variational HF calculation of the SF, for two systems. 
The SF is essentially correct beyond a system-size de- 
pendent minimum kg. At fc < fco, the SF levels off and 
at fc = approaches a value that equals the error in the 
xc hole. The crossover ko obviously goes to zero as the 
system size increases. Since our xc-hole is expected to 
be accurate inside a sphere of radius uq, one expects the 
SF to be accurate beyond 27r/(2uo), 2uo being the char- 
acteristic length scale of the simulation cell. Indeed, a 
cutoff fco ~ tt/uq seems plausible. We have looked at 
the fco values at which our HF structure factor and the 
exact one start to diverge markedly, and we have found 
fco ~ 4.2/uo. This rough estimate implies an accurate 
HF xc hole within a radius ~ 3mo/4 from a given elec- 
tron. Using fco = 4.2/mo and letting our calculated HF 
structure factor go to linearly at smaller fc produces a 
new corrected estimate of Uxc shown in Fig. |31 There is 
no longer a systematic error due to the system size and 
the corrected SF exhibits fluctuations that are consider- 
ably smaller than those of the uncorrected SF. Statistical 
noise and a shell structure remain, of course. 

Figs. 14151 repeat the analysis for the correlated Slater- 
Jastrow wavefunction. Here we use the VMC SF of the 
largest system as the reference Interestingly, the 

VMC structure factor seems to remain correct at unex- 
pectedly small values of fc. We have found fco ^ 2.1/uo, 
corresponding to a correct xc hole up to a surprisingly 
large u = 1.5mo! Possibly, 5*0 = poses such a strong 
constraint on the relatively shapeless ^3 SF that the SF 
has little choice but to be accurate at "too small" fc espe- 
cially as in contrast to the HF case, the true interacting 
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SF is also quadratic for fc ^ The behavior of 

the SF around fc = is therefore quahtatively correct. 
Nevertheless, a correction is needed. We chose a simple 
scheme, multiplying S{k) by k/ko when k < k^. We did 
not use the RPA for interpolation as its region of validity 
seems to begin at values of k smaller than our smallest 
/co (see also Fig. 0J). By looking at the final convergence 
of the xc interaction energy Uxc with system size, we see 
again that the estimates for Uxc are noisy but essentially 
flat, which is a signature of Coulomb finite-size errors 
having been eliminated. 

QMC being inherently statistical in nature, we close 
with a discussion of error bars. Errors of the SF at dif- 
ferent k are correlated and so the direct evaluation of 
an error for an xc interaction energy Uxc derived in the 
way described here is non-trivial. However, observe that 
in a finite system the SF-based Uxc and the usual Uxc 
coincide. Even for a infinite system, as the simulation 
cell increases their values become more and more simi- 
lar as they contain similar if not identical information. 
Hence, it seems reasonable to assume that the error bars 
of the standard Uxc can be used for the SF data. We 
have evaluated both the SF-based and the standard Uxc 
for a homogeneous electron gas with 2 * 51 electrons, for 
10 statistically independent yet identical runs averaging 
over 1000 QMC steps. The estimates for Uxc differ by an 
offset, due to the different finite-size errors, nevertheless, 
the estimated standard deviations are similar: 0.00182 
using the Ewald interaction, 0.00216 using the uncor- 
rected SF, and 0.00173 using the corrected SF. These re- 
sults are consistent with an error bar of 0.0015 derived by 
blocking. The estimates for Uxc were correlated, with the 
correlation between the Ewald calculation and the un- 
corrected and corrected SF-based calculations being 0.64 
and 0.74 respectively. The two SF-based calculations had 



a correlation coefficient of 0.88. Thus the error-bars for 
the standard Ewald data can be used also for the SF 
based data. 

In conclusion, we have devised a new method to eval- 
uate QMC xc interaction energies Uxc that do not suf- 
fer from spurious interactions of electrons with per iodic 
copies of their xc hole. The method is robust[l3 and 
easy to implement. Applying our method to Slater and 
Slater- Jastrow type many-body wavefunctions of a ho- 
mogeneous electron gas, we have shown how to efficiently 
handle and eliminate residual Coulomb finite-size errors. 
Our approach is equally applicable to arbitrary inho- 
mogeneous many-electron systems. Spherical averaging 
reduces the information contained in any QMC system 
to a smooth one-dimensional curve. Each value of S{k) 
therefore contains more information (hence less statisti- 
cal noise) than S{k, —k) and no spherical self-averaging 
need be assumed. In future, we aim to apply our method 
to the case of Jellium surfaces and real solids. In the case 
of the homogeneous electron gas that we have considered 
here, we have shown that SF-based calculations of the 
xc interaction energy can be improved considerably by 
simply letting the SF go to the correct long-wavelength 
limit at fc = 0. For more complex inhomogeneous many- 
electron systems it might be advantageous, however, to 
splice together structure factors obtained from RPA at 
small k and QMC at larger k. 
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